Mud pulse telemetry systems and methods using receive array processing

ABSTRACT

Mud pulse telemetry systems and methods may employ a downhole pulser to encode a digital telemetry data stream as pressure fluctuations in a fluid flow stream. An arrangement of spatially separated sensors acquire pressure-responsive measurements at multiple positions within the plumbing of a drill rig. A receiver collects and digitizes the measurements from the spatially separated sensors and subjects them to a principal components analysis (PCA) to determine those one or more basis vectors associated the telemetry signal. The PCA process may employ decomposition of a spatial correlation/covariance matrix or a temporal-spatial correlation/ covariance tensor. The selected basis vector(s) are then used to obtain the telemetry signal for demodulation into the telemetry data stream.

CROSS-REFERENCE TO RELATED APPLICATIONS

This application claims priority to Provisional U.S. Application Ser.No. 61/684,731, titled “Mud Pulse Telemetry Systems and Methods UsingReceive Array Processing” and filed Aug. 18, 2012 by Victor J. Stolpman,which is hereby incorporated herein by reference.

BACKGROUND

In the exploration of oil and gas, hydrocarbon seeking entities oftenemploy Measurement While Drilling/Logging While Drilling (MWD/LWD)services to measure and collect drilling and/or formation information inreal time while drilling (as the name suggests). In their quest, theseentities seek out drilling service companies that specialize in suchservices as directional drilling, formation evaluation, and otheroil-well drilling/planning operations. These services companies willthen bring various electronic “down-hole tools” to a well site andinsert them into the drill string near the drilling bit to sensedrilling or rock properties, acquire, process, transmit and present thereal-time information as the well is being drilled.

In addition to providing MWD/LWD measurements, the downhole toolsfacilitate directional drilling. Directional drilling employs asurveying instrument that estimates the orientation and optionally theposition of the bit. The surveying instrument is coupled with a steeringmechanism that enables the driller to navigate the borehole relative tosubsurface regions of interest.

The assembled downhole tools are often collectively referred to as theBottom-Hole-Assembly (BHA) that includes a telemetry module forcommunicating with the surface. The BHA transmits a representation of aleast a subset of the measured/logged data to the surface, where it isprocessed and presented to a user as a log report. Often, this data isof a digital format and may subject to a variety of compressiontechniques. There are a variety of transmission methods including butnot limited to mud pulse telemetry, electro-magnetic telemetry, acoustictelemetry and wired-pipe telemetry. Of these, mud pulse telemetry hasproven most popular due to its reliability.

Mud pulse telemetry (MPT) exploits the drilling rig's plumbing system.In most drilling operations, a circulation pump circulates fluid througha drill string and out the drill bit into a borehole, where it returnsalong the annulus to the surface. This fluid (often called “mud” in theoilfield industry) may include water and/or oil and one of a pluralityof additional additives that may be inert or chemically reactive withother molecular compositions present within a borehole during drillingoperations. There are a multitude of motivations for pumping mud withone example being simply to remove earth materials from the borehole.

A service company may install at least one transducer/sensor within therig's plumbing system. The surface rig's plumbing system mechanicallyconnects the so-called mud-pump(s) with the drill-string which in turnscouples with a drill-bit within the borehole. The purpose of thesesurface transducers/sensors is to enable the acquisition of encodedwaveforms transmitted from a “pulser” down in the BHA. Thus, MPT systemsserve to communicate real-time information through a mud column within adrill string via a series of modulated pressure waves that are supposedto be detected by the surface transducers/sensors.

Unfortunately, signal energy dispersion in the fluid and noise sources(e.g., from the circulation pump) often hinder the operation of the MPTsystem. Various precautions may be employed such as usage of pumpdampeners (a.k.a. “de-surgers”), careful positioning of the transducer,and directional detection. Such precautions are usually helpful, yetthere remain opportunities to further enhance telemetry systemperformance.

BRIEF DESCRIPTION OF THE DRAWINGS

Accordingly, there are disclosed in the drawings and detaileddescription specific embodiments of mud pulse telemetry methods andsystems that employ receive array processing. In the drawings:

FIG. 1 shows an illustrative measurement-while-drilling (MWD)environment suitable for the disclosed systems and methods.

FIG. 2 shows various forms of positive and negative pulsers.

FIG. 3 is a flow diagram of an illustrative MPT method employing receivearray processing.

It should be understood, however, that the specific embodiments given inthe drawings and detailed description do not limit the disclosure. Onthe contrary, they provide the foundation for one of ordinary skill todiscern the alternative forms, equivalents, and modifications that areencompassed in the scope of the appended claims.

TERMINOLOGY

In the present application, the following acronyms may be employed. Theyare spelled out here for easy reference:

AR Auto Regressive

ED Eigenvalue Decomposition

MA Moving Average

PCA Principal Component Analysis

SVD Singular Value Decomposition

DETAILED DESCRIPTION

The following description relates to a variety of mud pulse telemetry(MPT) method, apparatus and system embodiments that enable MeasurementWhile Drilling (MWD) services with real-time data transfer from sensorsnear the bit within a borehole to a surface location. Although notlimited to surface systems, this disclosure does focus most of thediscussion on the receiver side configurations, but this does not implythat this disclosure is limited to surface systems. One skilled in theart will recognize non-surface system embodiments (e.g., downlinktelemetry systems) are readily derivable from the ensuing description.

In short, this disclosure describes MPT systems and methods that employarray processing techniques such as Principal Component Analysis (PCA)to a plurality of digitized measurements from a plurality of transducers(sensors) spatially separated along a drilling rig's plumbing. SingularValue Decomposition (SVD) and Eigenvalue Decomposition (ED) are examplesof ways to extract principal components that exploit 2nd orderstatistical similarities/differences (e.g. spatial and/or temporalcorrelations) between said digitized measurements of a physicalphenomenon (e.g. pressure, velocity, acceleration, shear, temperature,phase changes, and wall strain).

The sensors each provide a signal with an associated output signal rangethat defines the limits on a corresponding axis in a “signal space”. Theset of sensor measurements when taken together form coordinates of apoint in the signal space. The full collection of sensor measurementsets defines a subspace within the signal space. This subspace can beexpressed in a number of ways, but a particularly useful representationis the set of principal components, which can be loosely expressed asbeing the set of orthogonal directions that best align with thesubspace, in decreasing order of signal variance. Signal decompositionmethods such as ED and SVD attempt to identify these principalcomponents in the form of orthonormal basis vectors (pairwise orthogonalvectors with unit norm length). The unit norm definition may vary. Acommon norm definition is Euclidean norm, i.e. l2-norm. Other commonnorms may be the l1-norm and (for discrete measurements) the l∞-norm.Embodiments using SVD and ED for PCA processing will likely use theEuclidean norm.

When ordered in rank of importance, the principal components attempt tocapture the most variance possible within the signal subspace observedin said plurality of digitized measurements. When using SVD, this formof signal decomposition may sort the discovered basis vectors indescending order of “interest” using at least in part the magnitude ofeach basis direction's corresponding singular value. Some embodimentsmay want to reduce receiver complexity and choose to utilize a reducedset of features for detection purposes rather than process the entireset of observed digital measurements. These embodiments may discardand/or ignore the components contributing only a minor portion ofinformation. These so-called “minor-components” may correspond to basiscomponents related to singular value magnitudes below some predeterminedthreshold (e.g. small absolute aggregate value or a percentage of sum ofall singular values).

Some method embodiments may further account for temporal components ofthe measurements with spatial-temporal tensor(s) to extract directionsrelevant to both space and time simultaneously. This extraction is doneby raking across time and space to construct 1st and 2nd orderstatistics to serve as the basis for PCA.

FIG. 1 depicts an illustrative Mud Pulse Telemetry (MPT) apparatusembodiment and system embodiment in use at a typical drillinginstallation while operating a Measurement While Drilling (MWD) service.As illustrated, the typical drilling installation includes a drillingderrick 102 at the surface of the well. The derrick may be transportableand temporarily erected on location. The drilling derrick 102 supportsthe drill string 104 and BHA 106 via a hoist 108 and swivel 110. In theFIG. 1 example, the BHA 106 includes a pulser 112, a tool sensor 114,and a drill bit 116. The BHA may further include additional MWD tools,stabilizers, and/or drill collars or heavyweight drill pipe (HWDP) tostiffen the BHA and add additional weight to aid with keeping the drillbit “on-bottom”.

The hoist 108 lowers drill string 104 through the rotary table 118 intothe casing 120 and beyond into the open borehole 122 until the bit 116reaches the bottom. The rotary table 118 turns the drill string 104 andbit 116 to extend the borehole through earth formations. If desired, adownhole mud motor can be employed to rotate the bit at a different ratethan the drill string.

Circulation pumps 124 take drilling fluid (“mud”) from a retention pit125 and circulate it through a feed pipe 126 to swivel 110 where itflows downward through the drill string interior as indicated by arrow128. Once the fluid reaches the bit 116, it exits through ports near thecutting elements to entrain and transport rock cuttings upward along theannulus as indicated by arrow 130. The fluid transports the cuttingsinto the retention pit 125 via return pipe 132. As the drilling mudcirculates through the drill bit, the drilling fluids functionadditionally as a bit coolant and lubrication extending the lifespan ofthe bit. Ideally, the weight and hydraulic pressure of the drillingfluid flow balances with the formation pressure to minimize fluid lossto the formation while still preventing an uncontrolled release offormation gases and fluids into the borehole, i.e. a “blowout.”

Pumps 124 are normally piston-based, causing a significant degree ofpressure variation due to the action of the pistons and valves. Apulsation dampener 134 is positioned along the feed pipe 126 toattenuate the (relatively) high-frequency variation, typically with onlya moderate degree of success. Downstream of the pulsation dampener, FIG.1 shows multiple transducers 136 that respond to pressure variation offluid in the feed pipe 126. The transducers 136 can be directly coupledto the fluid to physically respond to pressure variations, or coupled toa tubular housing the fluid flow to measure dimensional changesresulting from pressure variation in the flow stream. The transducerprovides a measurable reference signal (e.g. voltage, current, phase,position, etc.) sensitive to the pressure or temporal derivativethereof, i.e. dP(t)/dt, with a response that is proportional to withinan understood distortion (e.g., scalar gain, constant phase shift,time-shift, finite precision, etc.).

A transducer interface 138 converts the transducer response into anelectrical signal suitable for digitization and processing by receiver140. Receiver 140 may be dedicated MPT receiver electronics or a generalpurpose computer with a data acquisition card and suitable software forprocessing the acquired transducer signal(s). Among other things,receiver 140 may include circuitry, firmware, or software that performsPCA analysis of the acquired signals and employs selected basisvector(s) to extract a telemetry signal. The receiver may beincorporated as part of a computer system or coupled to a computersystem that provides a graphical user interface (GUI) to enable a userto view and optionally interactively analyze the telemetry data. Thecomputer may de-multiplex the telemetry data stream to obtaintool-specific logs which can be subjected to individual processing toconstruct and display the logs in an image format.

To communicate with the surface, a downhole “pulser” induces pressurefluctuations in the flow stream 128. The pressure fluctuations propagateupstream as pressure waves 142 until they reach the transducers 136.Information can be encoded into the pressure waves via modulation suchas frequency modulation, phase modulation, pulse position modulation,and pulse width modulation. Other suitable modulation schemes alsoexist. The chosen modulation scheme preferably provides sufficientdetection signal-to-noise ratdespite the attenuation, dispersion, andnoise effects introduced into the flow stream 128.

As part of the BHA 106, the down-hole pulser 112 may be mechanicallyand/or electrically coupled with additional down-hole sensors 114 thatmeasure, calculate and/or sense various conditions within or near thebottom of the borehole being drilled. The BHA may have an electricalpower source and inter-communicating control buses that facilitate thetransfer of data between BHA components. Not limited to the following,the electrical power source may be batteries and/or generator-basedderiving power from the flow of fluids via turbine or like mechanism.Likewise, not limited to the following, said control bus lines may be ofa metallic, conductive material for use with electrical systems and/ordielectric material when used with optical sources. FIG. 1 illustrates asingle downhole tool sensor coupled with a pulser, but those skilled inthe art understand MWD BHA configurations may have a multitude of toolsabove and/or below a pulser and may utilize more than one communicationmedia, e.g. mud pulse and electromagnetic telemetry.

A downhole controller may be included in the BHA with electroniccircuitry that collects from the various sensors measured formationevaluation values such as density of rock formation, pressure of thedrilling fluid, and gamma ray readings, and resistivity of rockformation. Additional measurements may include directional informationsuch as but not limited to inclination, tool-face, azimuthal, and/orsurveys. The controller may include an encoding module (e.g., in theform of circuitry or a programmable processor executing software in anassociated memory device) that encodes the collected information as adata stream for transmission by the pulser 112.

The pulser 112 actuates a valve at least in part to encoding themeasurement data stream as pressure modulations of the flow stream. FIG.2A shows a first illustrative pulser implementation having a valve orvariable flow restrictor formed from a circular, fan-like stator 201having multiple fan blades/fins extending radially from a central hub,and a similarly shaped rotor 202 that can oscillate with respect to the(stationary) stator 201. In this implementation, the valve is said to beclosed when the relative alignment of the stator and rotor finsmaximally restricts fluid flow (by misaligning the openings betweenblades). It is said to be open when the relative alignment of the statorand rotor fins minimally restricts fluid flow (by aligning the openingsbetween blades).

The valve is coupled serially within the fluid column to restrict (whenclosed) or ease (when open) the flow of fluid through the valve towardsthe drill-bit. When the valve is closed, a pressure build up occurswithin the fluid on the source side creating a positive pressure changethat propagates up to the surface. A subsequent opening of the valveenables the upstream pressure to drop to its previous pressure. Thus asthe rotor 202 oscillates, the valve creates a periodic pressurepulsation that is amenable to frequency and phase modulation.

FIG. 2B shows a second illustrative pulser implementation having aspinning rotor 204 in place of an oscillating rotor 202. As before, thealternation between alignment and misalignment of the openings betweenthe blades produces a periodic pressure pulsation that can be frequencyand phase modulated. A spinning rotor may offer better frequencystability at the expense of a more limited modulation range.

FIG. 2C shows a third illustrative pulser implementation having a floworifice 206 and a poppet 208 that moves relative to the orifice torestrict (when closed) and ease (when opened) the flow of fluid throughthe valve. A closing and re-opening of the valve (also referred to as amomentary closing of the valve) generates an upgoing pressure pulse(“positive pulse”).

FIG. 2D shows a fourth illustrative pulser implementation, which isoften termed a “negative pulser”. This pulser configuration includes abypass valve to vent fluid from the drill string bore into the annulus,thereby bypassing the drill bit. This venting of drilling fluid producesa pressure drop (i.e. a negative pressure change) within the drillstring's fluid column. FIG. 2D shows a valve seat and gate 210configuration. The gate 210 moves relative to the seat to close thevalve (i.e., restrict fluid flow into the annulus) and open the valve(permit fluid flow into the annulus). After closing the valve, the fluidpressure immediately rises in the drill-string column towards thesteady-state pressure prior to the valve's opening of the valve. As thename suggests, this opening and closing actuation of the valve creates anegative pulse that propagates throughout the column of drilling fluid.

In the configurations of FIGS. 2C and 2D, the valve is controllable togenerate individual pressure pulses that propagate to the surface,enabling the use of pulse width modulation and pulse positionmodulation. The modulation (whether frequency, phase, pulse position,pulse width, or some other form of modulation) is handled by thereceiver after the pressure variation signals have been acquired via thetransducer(s) 136.

Conventional strain gauge sensors may serve as transducers 136 toprovide a measurable reference, e.g. 4-20 mA current, proportional tothe mechanical fluid pressure present at the coupling point, i.e. P(t),by being directly coupled to the drilling fluid flow. Alternatively,such strain gauges could be employed to measure the strains that therig's plumbing undergoes when a mud pulse is present and/or absent.Examples of manufacturers of said sensors include but not limited toHoneywell and Rosemount (Emerson Electric affiliated).

Alternatively, as shown in FIG. 1, transducers 136 may respond to thetemporal derivative of the pressure signal at each coupling point, i.e.,dP(t)/dt, or the derivative of a commensurate strain in the plumbingthat is proportional to the pressure signal derivative. For the formermeasurement, direct coupling of the transducer to the fluid flow can beused. For the latter, the transducer can be coupled to the surface of atubular in the drill rig's plumbing (e.g., feed pipe 126). In theillustrated embodiment, transducers 136 each include an optical fiberwinding on the feed pipe to measure the strains via small changes in thefeed pipe dimensions.

In a most general sense, there may be different sensor types atdifferent positions withing the rig plumbing. A strain-gauge pressuresensor may have a direct connection with the pressurized fluid via aplumbed fitting. Alternatively, fiber optic or strain gauge sensors maybe attached to the exterior of the plumbing to sense pressure-induceddeformations, or in some cases, the time derivative of suchdeformations. Ultrasonic flow sensors may be employed to measurepressure-induced velocity variations, flow shear variations, or evendensity variations. There may be other flow properties and associatedsensor types that measure flow properties indicative of MPT signals.

Moreover, at least some disclosed embodiments contemplate coupling oneor more sensors at each of a plurality of connection points withindrilling rig's plumbing (e.g. stand-pipe, pumps, manifold, etc.). Unlessmultiple sensors are at the same exact location, each sensor will see adifferent waveform with a different viewpoint of the pressure waves.This disparity is partially due to location with respect to dispersivecharacteristics of the drilling rig's plumbing, but also due to relativelocation to a noise source (e.g. pumps), and waveformfrequency/wavelength-dependent components that may interfereconstructively or destructively depending on relative location ofsignificant reflectors and the multi-path effects of such reflectors.

Often, at least one strain gauge sensors is placed on the so-called“standpipe”, which is typically the farthest feasible mounting pointfrom the pump noise sources. (The standpipe is a vertical pipe thatconnects to a flexible high-pressure hose to deliver the fluid flow tothe swivel 110.) Strain gauge sensors typically provide a 4-20 mAreference current that after passing through a resistive load may bedigitized as a sampled voltage or current measurement. This digitizedmeasurement is proportional to pressure one would see in the standpipe.Other types of sensors that offer great flexibility in location arenon-invasive fiber optic sensors that sense minute changes in pipecircumference via phase changes within the light beams traveling withinthem. Appropriately configured fiber optic sensors yield measurementsproportional to the time derivative of pressure within the pipe.

As suggested earlier here within, this circulating column of drillingfluid flowing through the drill string may function as a transmissionmedium for encoded pressure pulses modulated via a pulser valve'sposition and/or orientation. These pressure pulses carry MWD informationfrom the BHA to the surface (and sometimes from the surface to the BHAin the case of a “downlink”). In the depicted system embodiment, thepulser and downhole sensors (a.k.a. “downhole tools”) are part of a mudpulse telemetry system, and surface systems may monitor sensor outputssensitive to the temporal differential pressure changes at each couplingpoint of the sensors and a tubular housing, e.g. the drilling rig'splumbing. Such embodiments may be observe the pressure either directlyfrom the drilling fluid and/or via small pipe diameter changes.

Additionally, the BHA may have an electrical power source andinter-communicating control buses that facilitate the transfer of databetween tools and pulser. Not limited to the following, the electricalpower source may be batteries and/or generator-based deriving power fromthe flow of fluids via turbine or like mechanism. Likewise, not limitedto the following, said control bus lines may be of a metallic,conductive material for use with electrical systems and/or dielectricmaterial when used with optical sources. FIG. 1 illustrates a singledownhole tool sensor coupled with a pulser, but those skilled in the artunderstand MWD BHA configurations may have a multitude of tools aboveand/or below a pulser and may utilize more than one communication media,e.g. mud pulse and electromagnetic telemetry.

With this context, then, we turn now to a more detailed discussion ofthe proposed array processing techniques. Specifically, this disclosuredescribes and applies Principal Component Analysis (PCA) to a pluralityof digitized measurements from a plurality of transducers (sensors)spatially separated along plumbing that transports the fluid flow.

Singular Value Decomposition (SVD) and Eigenvalue Decomposition (ED) areexamples extracting principal components that exploit 2nd orderstatistical similarities/difference (e.g. spatial and/or temporalcorrelations) between digitized measurements produced by the sensors. Atleast some receiver embodiments observe the plurality of digitizedmeasurements associated with the pressure waves traveling in theplumbing at a well site and process said measurements to estimatevarious 1st and 2nd order statistics, (e.g. sample mean, samplevariance, sample standard deviation, Auto Regressive (AR) approximationof ergodic and/or spatial means and correlations, Moving Average (MA)sample means and sample covariances within a sliding window related atleast in part to the MA order). Other embodiments may approximate theprobability density functions (e.g. joint, conditional, marginal) andcompute the various expected 1st and 2^(nd) order statistics using aprocessor.

In at least some embodiments, software configures the processor toconstruct at least one array of statistical estimates within a firstmemory storage, e.g. a covariance matrix and/or an estimated correlationmatrix (this latter may assume the means are zero or have been made zerovia a high pass filter, a mean subtraction step, or other mean removaldevice). These embodiments may further determine at least one featurerelated to at least one principle component. Some embodiments mayfurther determine a set of principle component vectors that form a basisfor a range spanned by the digital measurements. That is, a processoroperating from a set of instructions may decompose thecovariance/correlation matrices into additional matrices and or vectors.These matrices (and vectors) may include at least one unitary matrixand/or at least one set of singular values stored in memory.

When ordered in rank of importance, the principal components attempt tocapture the most variance possible within the observed signal subspaceprovided by the digitized measurements. When using SVD, this form ofsignal decomposition may sort the discovered basis vectors in descendingorder of “interest” using at least in part the magnitude of each basisdirection's corresponding singular value. Similarly, other embodimentsmay follow a similar path of Eigenvalue Decomposition and sort a seteigenvector directions (also an orthonormal basis) according to theircorresponding eigenvalues.

Some receiver embodiments may attempt to reduce processing complexity byemploying a reduced set of sensors for detection purposes rather thanprocess the entire set of observed digital measurements. Alternatively,or in addition, these embodiments may discard and/or ignore the thoseprincipal components that contribute only a minor portion ofinformation. These so-called “minor-components” may correspond to basiscomponents related to singular value magnitudes below some predeterminedthreshold (e.g. small absolute aggregate value or a percentage of sum ofall singular values).

Due to wave dispersion and the finite speed of wave propagation throughthe flow stream, temporal correlations are also expected between thesensor measurements. To identify and exploit such correlations, thereceiver “rakes” the signals relative to each other with varying timeshifts to construct 1st and 2nd order sample statistics which as beforeare subjected to PCA to extract principle components of thetime-extended signal subspace.

Before addressing the PCA analysis in great detail, consider FIG. 3which summarizes an illustrative MPT method. In block 302 the BHAobtains formation measurements and employs them to form a telemetrydatastream. The datastream may comply with a communications protocolthat, depending on the chosen implementation, includes compression,error correction coding, and framing with synchronization/trainingheaders and checksums. The communications protocol may further providechannel coding and modulation requirements, such as frequency shiftkeying, phase shift keying, on-off keying, amplitude modulation, pulsewidth modulation, pulse position modulation, differential encoding,channel pre-coding, and run-length limitations.

In block 304, the pulser modulates the flow stream to generate pressurefluctuations representing the telemetry datastream. The pressurefluctuations propagate along the flow string to the array of sensors.

In block 306, the surface receiver collects measurements from the arrayof sensors in the rig plumbing. As previously discussed, the sensors maybe spatially separated and in some cases may be of different types tomeasure different flow properties that reflect the pressurefluctuations. The measurements may be digitized and buffered in memory.In some embodiments, the measurements may be stored for laterprocessing.

In block 308, the receiver applies PCA to the measurements to obtain thebasis vectors for uncorrelated sources of variation in the measurementset. The measurement set may be expected to have basis vectorsassociated with the telemetry signal as well as pump noise, pipevibration modes, and other noise sources. The basis vectors can bedetermined by post-processing the measurement data, but it may bepreferred to perform the PCA in an iterative or adaptive fashion toavoid introducing unwanted delay into the telemetry process. As analternative approach, a training phase may be employed to initiallyidentify the basis vectors for subsequent use. Such training may beperformed during receiver setup or periodically throughout the usage ofthe receiver.

In block 310, the receiver identifies those one or more basis vectorsassociated with the telemetry signal, and employs those basis vectors toextract the telemetry signal from the signal set. The identification canbe performed in a number of ways, including pattern matching or evenjust trial and error. For example, telemetry signal may be expected tohave a recognizable pattern (particularly if the communications protocolprovides for a synchronization/training header), so examination of thosesignals that result from the application of each basis vector to the setof measurements generally indicates which basis vector(s) extracts thetelemetry signal. Also, the success of the subsequent method steps willbe contingent on the proper basis vector selection, so the basis vectoryielding the lowest error rate or highest estimated detectionsignal-to-noise ratmay be selected as the appropriate basis vector.

In block 312, the receiver equalizes, demodulates, and decodes theextracted telemetry signal to obtain the telemetry data stream. Theequalization compensates for residual signal dispersion, while thedemodulation and decoding steps reverse the channel coding andmodulation steps, and may further include de-framing, error correctionand decompression.

In block 314, the system (either receiver 140 or one of the downstreamprocessing systems) de-multiplexes the data stream to obtain and storedata for each of the tools. The extracted log data is then subjected totool-specific processing to construct logs for display to a user.

We now turn to a more detailed discussion of block 308. PrincipalComponent Analysis (PCA) is a classic statistical data analysistechnique. The main idea of PCA is a reduction of dimensionality of agiven multi-variate data set from M to K. Using vector representation,let y=(y₀, y₁, . . . , y_(M−1))^(T) be a M×1 vector of M observationsincluding a set of at least K interrelated multi-variables, say u=(u₀,u₁, . . . , u_(K−1))^(T), where u is a K×1 vector and K≦M. Often therandom elements of u cannot be directly observed unobstructed or itsdimensionality K even be known a priori. Nevertheless, the goal is thento estimate this smaller set of “features,” say û of dimensionality K,from the possibly noise corrupted sampling of y that may statisticallyrepresent the K most interesting variables within the M observations,i.e. the principal components with the K largest variances.

This so-called reduction to (or extraction of) “principal” features fromwithin a stochastic process may employ one or more orthogonal (ororthogonal-like) transformations of y to reduce the multi-variable set'sdimensionality, i.e. u=Gy where y∈

^(M×1) and G∈

^(k×M) whose M columns are pairwise orthogonal. For example, using thewell-known dot-product as the inner product, the ith and jth columnvectors of G, g_(i) and g_(j) respectfully, are said to be pairwiseorthogonal if g_(i) ^(T)g_(j)=0 if i<>j and g_(i) ^(T)g_(j)=1 if andonly if i=j, ∀i. If in addition K=M, then G is said to be an orthogonaltransform, and thus, G^(T)G=I_((M×M)) where I_((M×M)) is the M×Midentity matrix.

Similarly, elements of y may include complex numbers taken from Mdimensional complex number space, i.e. y∈

^(M×1) and the orthogonal transform condition becomes G^(H)G=I_((M×M))where G^(H) is the conjugate-transpose of G, sometimes called theHermitian Transpose where G∈

^(M×M).

Throughout the remaining disclosure, the author will continue describingthe various embodiments of PCA assuming real number elements coming from

-space and will only address complex number space,

-space, when appropriate. It is sufficient to say that these conceptseasily extend to observations with complex number elements.

Mathematically, SVD is the matrix factorization of any given (M×N)matrix A. In the event of 1-dimensional real value measurements, we canwrite the SVD matrix factorization as

A=ΣV^(T)

where U is the (M×M) left unitary matrix factor, V^(T) is a (N×N) leftunitary matrix factor, and Σ is a diagonal matrix of (M×N) containingsingular values on its main diagonal, i.e.

UU ^(T) =U ^(TU) =I _((M×M))

VV ^(T) =V ^(T) V=I _((N×N))

Σ=diag(σ₁, σ₂, . . . , σ_(min(M,N)), 0, . . . , 0)

In embodiments that acquire a 2-dimensional sensor reading (e.g. phaseand amplitude, I and Q plane), one can then use the complex value formof SVD, e.g.

A=UΣV^(H)

Mathematically, eigenvalue decomposition (ED) is the matrixfactorization of a symmetric positive-semi-definite (p.s.d.) (M×M)matrix A. In the event of 1-dimensional real value measurements, we canwrite the ED matrix factorization as:

A=UΛU^(T)

where U is the (M×M) unitary matrix factor having eigenvectors and Λ isa diagonal matrix of (M×M) containing eigenvalues on its main diagonal,i.e.

UU ^(T) =U ^(T) U=I _((M×M))

Λ=diag(λ₁, λ₂, . . . , λ_(M))

In embodiments that acquire a 2-dimensional sensor reading (e.g. phaseand amplitude, I and Q plane), one can then use the complex value formof ED, e.g.

A=UΛU^(H)

ED is essentially a special case of SVD suitable for factoring squarep.s.d. matrices. Due to symmetry, ED's complexity is often lower thanSVD's complexity when M=N, but due to numerical issues, solutions forSVD are often more stable (i.e. less chance of finite precision issues)and thus used more often.

Turning now to measurements acquired by an array of spaced-apartsensors, let x[n]=(x[n], x[n−1], . . . , x[n−N+1])^(T) include N realvalues, i.e. x[n] is a (N×1) sample vector of the transmitted pressurewave signal as a function of time n. Using a general linear channelmodel, we can write for real valued output (likewise a similarconstruction exists for complex valued signals):

Y[n]=(y₁ [n], y ₂ [n], . . . , y _(NRX) [n]) ^(T)

y[n]=Hx[n]+η[n]

y[n]=UΣV ^(T) x[n]+η[n]

where y_(i)[n] is the ith received sensor measurement at time n, N_(RX)is the number of sensors, and H=UΣV^(T) is our channel's linear modeldecomposed.

By left multiplying by an estimate of U^(T), the processor will attemptto reverse a portion of the effects of the channel,

U ^(T) y[n]=U ^(T) Hx[n]+U ^(T) η[n]

U ^(T) y[n]=U ^(T) UΣV ^(T) x[n]+{tilde over (η)}[n]

U ^(T) y[n]=ΣV ^(T) x[n]+{tilde over (η)}[n]

The left multiply by U rotates the spatial signal space via a unitarytransform, i.e. invertible transform so no information loss occurs.Various embodiments estimate U via the SVD and/or ED of a discreteestimate of R_(yy)(τ_(k)), namely {circumflex over (R)}_(y,t)[τ_(k)]using either a Moving Average (MA) or Auto Regressive (AR) process. Theright unitary matrix factor V^(T) can then be readily derived from ameasure of the channel response H and used as the basis for extractingtransmitted signals.

To estimate R_(yy)(τ_(k)) and μ_(y), write the continuous and sampleddigitized measurements of N_(RX) receivers in vector form

y(t)=(y ₁(t), y ₂(t), . . . , y _(NRX)(t))^(T)

y[n]=(y ₁ [n], y ₂ [n], . . . , y _(NRX) [n])^(T)

The expected mean of the N_(RX) received signals can be written invector form as

μ_(yi) =E{y _(i)(t)}

μ_(y) =E{y(t)}

μ_(y) =E{(y ₁(t), y ₂(t), . . . , y _(NRX)(t))^(T)}

μ_(y)=(E{y ₁(t)}, E{y ₂(t)}, . . . , E{y _(NRX)(t)})^(T)

μ_(y)=(μ_(y1), μ_(y2), . . . , μ_(yNRX))^(T)

along with the covariance matrix

R_(yi, yj)(τ_(k)) = E{(y_(i)(t) − μ_(yi))(y_(j)(t − τ_(k)) − μ_(yj))^(T)}${R_{y,y}\left( \tau_{k} \right)} = {{E\left\{ {\left( {{y(t)} - \mu_{y}} \right)\left( {{y\left( {t - \tau_{k}} \right)} - \mu_{y}} \right)^{T}} \right\}} = \begin{pmatrix}{R_{{y\; 1},{y\; 1}}\left( \tau_{k} \right)} & {R_{{y\; 1},{y\; 2}}\left( \tau_{k} \right)} & \ldots & {R_{{y\; 1},y_{N_{RX}}}\left( \tau_{k} \right)} \\{R_{{y\; 2},{y\; 1}}\left( \tau_{k} \right)} & {R_{{y\; 2},{y\; 2}}\left( \tau_{k} \right)} & \ldots & {R_{{y\; 2},y_{N_{RX}}}\left( \tau_{k} \right)} \\\vdots & \vdots & \ddots & \vdots \\{R_{y_{N_{RX}},{y\; 1}}\left( \tau_{k} \right)} & {R_{y_{N_{RX}},{y\; 2}}\left( \tau_{k} \right)} & \ldots & {R_{y_{N_{RX}},y_{N_{RX}}}\left( \tau_{k} \right)}\end{pmatrix}}$

Assuming the means are zero or made to be zero by various embodiments,we can write:

R _(y,y)(τ_(k))=UΣ ² U ^(T) +R _(η,η)(τ_(k)); and

R _(y,y)(τ_(k))=UΣ ² U ^(T)+σ² Iδ(τ_(k));

if the noise is white. Thus using an estimate of R_(y,y)(τ_(k)) mayprovide an estimate of U and hence U^(T)y[n].

Moving Average or Auto Regressive approaches may be preferred for samplestatistics estimation, but other approaches could also be used.

Averaging over a sliding window of time and computing the so-calledMoving Average (MA) estimate is one methods for estimating the various1st and 2nd order statistics and constructing the mean vectors andcovariance/correlation matrices. The MA sample mean vector may bedefined as:

μ̂_(y)^(MA) = (μ̂_(y₁)^(MA), μ̂_(y₂)^(MA), …  , μ̂_(y_(N_(RX)))^(MA))^(T),

where the elements are the estimated means of the ith sensor:

${{\hat{\mu}}_{yi}^{MA} = {{\frac{1}{N_{MA}}{\sum\limits_{n = 0}^{N_{MA} - 1}\; {{y_{i}\lbrack n\rbrack}\mspace{14mu} {where}\mspace{14mu} i}}} = 1}},2,\ldots \mspace{14mu},N_{RX}$

Similarly, the MA sample covariance for τ_(k) may be defined as:

${{\hat{R}}_{y,y}^{MA}\left\lbrack \tau_{k} \right\rbrack} = \begin{pmatrix}{{\hat{R}}_{{y\; 1},{y\; 1}}^{MA}\left\lbrack {n,\tau_{k}} \right\rbrack} & {{\hat{R}}_{{y\; 1},{y\; 2}}^{MA}\left\lbrack {n,\tau_{k}} \right\rbrack} & \ldots & {{\hat{R}}_{{y\; 1},y_{N_{RX}}}^{MA}\left\lbrack {n,\tau_{k}} \right\rbrack} \\{{\hat{R}}_{{y\; 2},{y\; 1}}^{MA}\left\lbrack {n,\tau_{k}} \right\rbrack} & {{\hat{R}}_{{y\; 2},{y\; 2}}^{MA}\left\lbrack {n,\tau_{k}} \right\rbrack} & \ldots & {{\hat{R}}_{{y\; 2},y_{N_{RX}}}^{MA}\left\lbrack {n,\tau_{k}} \right\rbrack} \\\vdots & \vdots & \ddots & \vdots \\{{\hat{R}}_{y_{N_{RX}},{y\; 1}}^{MA}\left\lbrack {n,\tau_{k}} \right\rbrack} & {{\hat{R}}_{y_{N_{RX}},{y\; 2}}^{MA}\left\lbrack {n,\tau_{k}} \right\rbrack} & \ldots & {{\hat{R}}_{y_{N_{RX}},y_{N_{RX}}}^{MA}\left\lbrack {n,\tau_{k}} \right\rbrack}\end{pmatrix}$

and have its elements estimated linearly by

${{\hat{R}}_{{yi},{yj}}^{MA}\left\lbrack \tau_{k} \right\rbrack} = {\frac{1}{N_{MA}}{\sum\limits_{n = 0}^{N_{MA} - 1}\; {\left( {{y_{i}\lbrack n\rbrack} - {\hat{\mu}}_{yi}^{MA}} \right)\left( {{y_{j}\left\lbrack {n - \tau_{k}} \right\rbrack} - {\hat{\mu}}_{yj}^{MA}} \right)^{*}}}}$where  i = 1, 2, …  , N_(RX)

where τ_(k) here is an integer relating a sample separation betweensensor i and sensor j measurements.

In essence, Auto Regressive (AR) Estimates average over an exponentiallydecaying window as a function of time to estimate the various 1st and2nd order statistics and constructing the mean vectors andcovariance/correlation matrices. The AR sample mean vector may bedefined as:

{circumflex over (μ)}_(y) ₁ ^(AR) [n]=α{circumflex over (μ)} _(y) ₁^(AR) [n−1]+(1−α)y _(i) [n] where 0<α<1

where the elements are

$\begin{matrix}{{{\hat{\mu}}_{y}^{AR}\lbrack n\rbrack} = {{{\alpha {{\hat{\mu}}_{y}^{AR}\left\lbrack {n - 1} \right\rbrack}} + {\left( {1 - \alpha} \right){y\lbrack n\rbrack}\mspace{14mu} {where}\mspace{14mu} 0}} < \alpha < 1}} \\{= \left( {{{\hat{\mu}}_{y\; 1}^{AR}\lbrack n\rbrack},{{\hat{\mu}}_{y\; 2}^{AR}\lbrack n\rbrack},\ldots \mspace{14mu},{{\hat{\mu}}_{y_{N_{RX}}}^{AR}\lbrack n\rbrack}} \right)^{T}}\end{matrix}$

where n is the current integer time sample and 0<α<1. The AR samplecovariance for τ_(k) at time sample n may be defined as:

${{\hat{R}}_{y,y}^{AR}\left\lbrack \tau_{k} \right\rbrack} = \begin{pmatrix}{{\hat{R}}_{{y\; 1},{y\; 1}}^{AR}\left\lbrack {n,\tau_{k}} \right\rbrack} & {{\hat{R}}_{{y\; 1},{y\; 2}}^{AR}\left\lbrack {n,\tau_{k}} \right\rbrack} & \ldots & {{\hat{R}}_{{y\; 1},y_{N_{RX}}}^{AR}\left\lbrack {n,\tau_{k}} \right\rbrack} \\{{\hat{R}}_{{y\; 2},{y\; 1}}^{AR}\left\lbrack {n,\tau_{k}} \right\rbrack} & {{\hat{R}}_{{y\; 2},{y\; 2}}^{AR}\left\lbrack {n,\tau_{k}} \right\rbrack} & \ldots & {{\hat{R}}_{{y\; 2},y_{N_{RX}}}^{AR}\left\lbrack {n,\tau_{k}} \right\rbrack} \\\vdots & \vdots & \ddots & \vdots \\{{\hat{R}}_{y_{N_{RX}},{y\; 1}}^{AR}\left\lbrack {n,\tau_{k}} \right\rbrack} & {{\hat{R}}_{y_{N_{RX}},{y\; 2}}^{AR}\left\lbrack {n,\tau_{k}} \right\rbrack} & \ldots & {{\hat{R}}_{y_{N_{RX}},y_{N_{RX}}}^{AR}\left\lbrack {n,\tau_{k}} \right\rbrack}\end{pmatrix}$

having its elements estimated recursively by

{circumflex over (R)} _(yi,yj) ^(AR) [n,τ _(k) ]=α{circumflex over (R)}_(yi,yj) ^(AR) [n−1,τ_(k)]+(1−α)(y _(i) [n]−{circumflex over (μ)} _(yi)^(AR))(y _(j) [n−τ _(k)]−{circumflex over (μ)}_(yj) ^(AR))* where i= 1,2, . . ., N _(RX)

where τ_(k) here is an integer relating a sample count separationbetween sensor i and sensor j measurements, n is the current integertime sample and 0<α<1.

PCA can be performed with a Temporal-Spatial Tensor Array. Using similarnotations, we may describe the received continuous and discrete timetemporal-spatial tensor as:

y(t,{right arrow over (τ)})=(y ₁(t), . . . , y _(N) _(RX) (t), y ₁(t−τ₁), . . . , y _(N) _(RX) (t−τ ₁), . . . , y ₁(t−τ _(K)), . . . , y _(N)_(RX) (t−τ _(K)))^(T)

y[t,{right arrow over (τ)}]=(y ₁ [n], . . . , y _(N) _(RX) [n], y ₁ [n−τ₁ ], . . . , y _(N) _(RX) [t−τ ₁ ], . . . . , y _(i) [n−τ _(K) ], . . ., y _(N) _(RX) [n−τ _(K)]))^(T)

where {right arrow over (τ)}=(τ₁, τ₂, . . . , τ_(NRX)). Like thespatial-only case, these temporal-spatial tensor embodiments mayestimate the matrix using discrete sampled statistic estimates ofR_(yy)({right arrow over (τ)}) and μ_(y)({right arrow over (τ)}) namely{circumflex over (R)}_(yy)[{right arrow over (τ)}] and {circumflex over(μ)}_(y)[{right arrow over (τ)}], respectively.

As with the spatial-only case, various receiver embodiments may estimatethe following mean and vector

μ_(y)({right arrow over (τ)})=E{y(t,{right arrow over (τ)})}

R _(yy)({right arrow over (τ)})=E{(y(t,{right arrow over(τ)})−μ_(y)({right arrow over (τ)}))(y(t,{right arrow over(τ)})=μ_(y)({right arrow over (τ)}))^(T)}

Via MA and AR models in a similar fashion, and estimate a similar U.Thus using an estimate of R_(yy)({right arrow over (τ)}) variousembodiments may estimate U from R_(y,y) ^(MA) ({right arrow over (τ)})and/or R_(y,y) ^(AR)({right arrow over (τ)}) for a temporal-spatialtensor unitary transform and various embodiments may then use U^(T)y[n,{right arrow over (τ)}] to extract the principle component features anddata from both time and spatially distributed measurements.

Embodiments employing SVD and/or ED may determine a first unitary matrixhaving orthonormal basis vectors from a plurality of possible unitarymatrices (a set) with each matrix in the set having sign permuted basisvectors of the first said determined unitary matrix. In the event thatsingular values (and/or eigenvalues) corresponding to two or more basisvectors (and/or eigenvectors) are identical, the plurality set mayfurther have unitary matrices that are permutations of columns (and/orrows) unitary matrices in addition to the sign permutation mentionedearlier.

In practice, the drilling conditions and channel conditions may varytemporally. In order to adapt to the new conditions withoutsimultaneously suppressing polarity changes attributable to transmitteddata or introducing extraneous polarity changes, at least someembodiments would track any changes in said estimated statistics derivedfrom samplings of said plurality of digitized measurements using anynumber of estimation techniques (e.g. MA and AR); construct a new arrayof updated estimated statistics; determine principal components withinsaid plurality of digitized measurements; and detect/decode messagesfrom a pulsing device. These embodiments may further include decomposingsaid new array of updated estimated statistics using at least in partSVD and/or ED. The basis vectors from consecutive decompositions may becompared and polarity changes applied to the later decomposition toachieve minimized deviation angles between corresponding basis vectors.

Realtime decomposing embodiments may process the current vector sampleof digitized measurements and/or any vector sample of digitizedmeasurements thereafter. Additional embodiments may further includebuffering/storing any historic digitized measurements and then arrayprocess said historic measurements with a decomposition occurring afterthe historic data in time (“post-processing” or “replay”). Alternativelyto determine statistical changes, embodiments may simply decomposeupdated estimated statistics on a periodic basis. Nevertheless, ineither case there may still exist the sign and order ambiguity due tothe aforementioned permutations, and the waveform and order of outputmay permute and/or sign change with each new array decomposition ofupdated estimated statistics.

In order to address this issue, some embodiments may further include aGraphical User Interface (GUI) wherein said GUI may include a “LOCK”and/or “UNLOCK” control feature. This feature may also be implemented asa physical switch coupled to said processor wherein the physical switchand/or GUI control feature may stop channel statistic tracking and“LOCK” the last determined unitary transform operator for continued usein transforming said digitized measurements after an actuation event ofsaid switch and/or GUI control feature. At least some receiverembodiments may further include an additional control feature for usewith graphical displays, a “SIGN” flip GUI control for each transformeddigitized measurement waveform displayed on a monitor.

Other embodiments may address the sign flipping and order permutingissue by tracking the pairwise angle between each basis directionpairing of a current unitary matrix decomposition and at least oneprevious unitary matrix. These embodiments pair current and former basisvectors according to pairwise minimum angles and then sign adjust thebasis orientation (positive or negative of a basis vector stilldescribes the same subspace or axis). The complexity of this sign andorder tracking grows combinatorially until some basis vector pairingsmay no longer be practically tracked. In these cases, embodiments maydrop the sign and order tracking for the minor components.

One illustrative system employing a surface detection system with asingle pressure transducer and a matched filter (without PCA) to receivesignals from a negative pulser pulsing at or near 9000 ft measureddepth, provided a 62% sample probability of correctly receiving data.Under identical in situ conditions this illustrative system using asurface detection system with four pressure transducers distributeduniformly along the drill rig plumbing between the pumps and thestandpipe and spatial PCA-based signal extraction followed by a matchedfilter, provided a 97% sample probability of correctly receiving data.Calculating data throughput as the product of data rate and probabilityof correct reception, the PCA-based detection system yielded a 56%improvement in data througput. Though not quite as positive, otherrelated experiments also yielded favorable results.

Numerous modifications, equivalents, and alternatives will becomeapparent to those skilled in the art once the above disclosure is fullyappreciated. For example, the foregoing description focuses on uplinkcommunication from the BHA to the surface, but this disclosure alsoapplies to downlink communication from the surface to the BHA. Suchdownlink communications may be used to convey commands and configurationparameters to control tool operations and/or steer the drill string.Some system embodiments may include sensors coupled to noise sources(such as the circulation pumps) to facilitate identification and removalof noise components from the digitized measurements. Examples of suchsensors may include vibration sensors, accelerometers, and strokeposition monitors. It is intended that the following claims beinterpreted to embrace all such modifications, equivalents, andalternatives where applicable.

What is claimed is:
 1. A mud pulse telemetry system comprising: a pulserthat transmits a digital data stream as pressure modulations of a flowstream; a plurality of spatially separated sensors each responding topressure modulations of the flow stream; and a receiver that processesdigitized signals from the plurality of sensors to determine a principalcomponent basis vector associated with said pressure modulations fromthe pulser, wherein the receiver further employs said principalcomponent basis vector to extract the digital data stream.
 2. The systemof claim 1, wherein the pulser employs one of differential pulseposition modulation, phase shift keying, and frequency shift keying, totransmit the digital data stream.
 3. The system of claim 1, wherein aspart of determining the principal component basis vector, the receiverperforms a decomposition of a matrix formed from sampled statisticsincluding at least one of mean, correlation, and covariance.
 4. Thesystem of claim 3, wherein the decomposition is an eigenvaluedecomposition or a singular value decomposition (SVD).
 5. The system ofclaim 1, wherein as part of determining the principal component basisvector, the receiver performs a decomposition of a temporal-spatialcovariance or correlation tensor.
 6. The system of claim 1, wherein aspart of determining the principal component basis vector, the receiverreduces a mean of the digitized signals.
 7. The system of claim 6,wherein as part of determining the principal component basis vector, thereceiver normalizes the digitized signals using respective estimatedstandard deviations.
 8. The system of claim 1, further comprising acomputer that de-multiplexes the data stream into tool-specific datasets for use in constructing and displaying logs to a user.
 9. Thesystem of claim 8, wherein the computer provides a graphical userinterface for displaying the logs to a user.
 10. The system of claim 9,wherein the graphical user interface enables the user to monitor thedecomposition and control polarity of each principal component basisvector.
 11. A mud pulse telemetry method that comprises: acquiringspatially separated measurements that are responsive to pressurefluctuations of a flow stream in a pipe; digitizing the measurements;processing the measurements to determine a principal component basisvector associated with a telemetry signal; and using the principalcomponent basis vector to extract a digital data stream from themeasurements.
 12. The method of claim 11, wherein said processingincludes: determining a matrix of sample statistic values; anddetermining an eigenvalue decomposition or singular value decompositionof the matrix.
 13. The method of claim 12, wherein said matrix is aspatial covariance or spatial correlation matrix.
 14. The method ofclaim 12, wherein said matrix is a temporal-spatial covariance orcorrelation tensor.
 15. The method of claim 12, wherein said processingfurther includes reducing means of the spatially-separated measurementsand normalizing standard deviations of the spatially-separatedmeasurements as part of determining the matrix.
 16. The method of claim12, further comprising repeating said processing to account for changesin the flow stream.
 17. The method of claim 16, further comprisingmonitoring relative polarities between principal component basis vectorsobtained with each repetition of the processing step and inverting basisvector polarities as needed to maintain consistent polarity of at leastsome basis vectors.
 18. The method of claim 11, wherein said extractinga digital data stream includes applying the basis vector to themeasurements to obtain a telemetry signal, and further includesdemodulating the telemetry signal to obtain the digital data stream. 19.The method of claim 18, wherein the telemetry signal represents thedigital data stream with one of pulse position modulation, phase shiftkeying, and frequency shift keying.
 20. The method of claim 11, furthercomprising deriving one or more borehole logs from the digital datastream and displaying said logs.